#lang racket


(require "math.rkt"
         "misc.rkt")


(provide (struct-out point)
         xyz
         x
         y
         z
         
         sph
         sph-r
         sph-th
         sph-fi
         
         cyl
         cyl-r
         cyl-th
         
         pol
         pol-r
         pol-th
         
         zero
         one
         ux
         uy
         uz
         
         p?
         zero?
         
         p=p
         p+p
         p-p
         p.p
         pxp
         p<\p
         
         p*n

         unit
         
         +x
         +y
         +z
         +xyz
         +sph
         +cyl
         +pol
         
         list<-point
         vector<-point)


(define (cartesian-x<-spherical r th fi)
  (* r (cos th) (sin fi)))

(define (cartesian-y<-spherical r th fi)
  (* r (sin th) (sin fi)))

(define (cartesian-z<-spherical r th fi)
  (* r (cos fi)))

(define (cartesian<-spherical r th fi)
  (list
   (cartesian-x<-spherical r th fi)
   (cartesian-y<-spherical r th fi)
   (cartesian-z<-spherical r th fi)))


(define (cartesian-x<-cylindrical r th z)
  (* r (cos th)))

(define (cartesian-y<-cylindrical r th z)
  (* r (sin th)))

(define (cartesian-z<-cylindrical r th z)
  z)

(define (cartesian<-cylindrical r th z)
  (list
   (cartesian-x<-cylindrical r th z)
   (cartesian-y<-cylindrical r th z)
   (cartesian-z<-cylindrical r th z)))


(define (point-print p port mode)
  (write-string "p" port)
  (write-string (format "~a" (list<-point p)) port))

(struct point (x y z)
  #:property prop:custom-write point-print)

(define (xyz x y z)
  (point x y z))

(define x point-x)
(define y point-y)
(define z point-z)


(define (sph r th fi)
  (apply xyz (cartesian<-spherical r th fi)))

(define (sph-r p)
  (sqrt (+ (^2 (x p)) (^2 (y p)) (^2 (z p)))))

(define (sph-th p)
  (if (= (x p) (y p) 0)
      0
      (atan (y p) (x p))))

(define (sph-fi p)
  (if (= (sph-r p) 0)
      (sph-r p)
      (acos (/ (z p) (sph-r p)))))


(define (cyl r th z)
  (apply xyz (cartesian<-cylindrical r th z)))

(define (cyl-r p)
  (sqrt (+ (^2 (x p)) (^2 (y p)))))

(define (cyl-th p)
  (sph-th p))

(define (cyl-z p)
  (z p))


(define (pol r th)
  (cyl r th 0))

(define (pol-r p)
  (cyl-r (struct-copy point p (z 0))))

(define (pol-th p)
  (cyl-th (struct-copy point p (z 0))))


(define zero (xyz 0 0 0))
(define one (xyz 1 1 1))
(define ux (xyz 1 0 0))
(define uy (xyz 0 1 0))
(define uz (xyz 0 0 1))

(define (p? p)
  (point? p))

(define (zero? p)
  (and
   (p? p)
   (= (x p) 0)
   (= (y p) 0)
   (= (z p) 0)))

(define (p=p p1 p2)
  (and
   (= (x p1) (x p2))
   (= (y p1) (y p2))
   (= (z p1) (z p2))))

(define (p+p p1 p2)
  (xyz
   (+ (x p1) (x p2))
   (+ (y p1) (y p2))
   (+ (z p1) (z p2))))

(define (p-p p1 p2)
  (xyz
   (- (x p1) (x p2))
   (- (y p1) (y p2))
   (- (z p1) (z p2))))

(define (p.p p1 p2)
  (+
   (* (x p1) (x p2))
   (* (y p1) (y p2))
   (* (z p1) (z p2))))

(define (pxp p1 p2)
  (xyz
   (-
    (* (y p1) (z p2))
    (* (z p1) (y p2)))
   (-
    (* (z p1) (x p2))
    (* (x p1) (z p2)))
   (-
    (* (x p1) (y p2))
    (* (y p1) (x p2)))))

(define (p<\p p1 p2 (normal #f))
  (define (aux)
    (acos
       (p.p
        (p*n p1 (/ (sph-r p1)))
        (p*n p2 (/ (sph-r p2))))))
  
  (define (continuous-aux)
    (let ((p (pxp p1 p2)))
      (if (or
           (p=p zero p)
           (> (p.p normal p) 0))
          (aux)
          (- 2pi (aux)))))
  
  (if (eq? normal #f)
      (aux)
      (continuous-aux)))


(define (p*n p n)
  (xyz
   (* (x p) n)
   (* (y p) n)
   (* (z p) n)))


(define (unit p)
  (p*n p (/ (sph-r p))))


(define (λxyz p (λp λid) (λx λid) (λy λid) (λz λid))
  (λp
   (λx (x p))
   (λy (y p))
   (λz (z p))))

(define (+x p n)
  (λxyz p xyz (cλ + n) λid λid))

(define (+y p n)
  (λxyz p xyz λid (cλ + n) λid))

(define (+z p n)
  (λxyz p xyz λid λid (cλ + n)))

(define (+xyz p n0 n1 n2)
  (λxyz p xyz (cλ + n0) (cλ + n1) (cλ + n2)))

(define (+sph p r th fi)
  (apply +xyz p (cartesian<-spherical r th fi)))

(define (+cyl p r th z)
  (apply +xyz p (cartesian<-cylindrical r th z)))

(define (+pol p r th)
  (+cyl p r th 0))


(define (list<-point p)
  (λxyz p list))

(define (vector<-point p)
  (λxyz p vector))